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Summary 


We  have  successfully  discriminated  the  quarry  blasts  and  earthquakes  at  local 
distance  from  the  observed  source  power  spectra  by  using  the  coda  Q  method.  To  extend 
this  discrimination  method  to  the  events  at  regional  distances,  we  have  developed  the 
energy  transfer  theory  for  a  mixture  of  Rayleigh  and  S  wave  propagation  in  a  random 
scattering  and  attenuation  medium.  Applying  this  theory  to  the  quarry  blast-earthquake 
problem  at  local  distance,  we  can  explain  the  observations  satisfactorily.  Our  study 
indicates  that  the  Rayleigh  waves'  scattering  become  important  for  the  events  with 
shallower  focal  depth  such  as  the  quarry  blasts  and  nuclear  explosions.  By  fitting  our 
theoretical  attenuation  curves  to  the  observed  ones  at  local  distances,  we  have  obtained 
the  S  and  Rayleigh  waves'  intrinsic  and  scattering  coefficients,  and  the  coupled  scattering 
coefficient  of  Rayleigh  with  S  waves.  Assuming  that  these  coefficients  are  valid  for  the 
events  at  regional  distances,  we  have  compared  the  predicted  and  observed  coda  envelope 
at  regional  distances  for  the  purpose  of  discrimination  between  the  earthquakes  and 
explosions  at  regional  distance.  Our  prediction  shows  that  the  Rayleigh  wave's 
contributions  become  more  significant  for  shallower  events  (explosions),  while  for  deeper 
events  (earthquakes)  the  Rayleigh  wave's  contributions  are  still  negligible.  Comparison 
between  predicted  and  the  data  of  CDSN  suggests  that  the  discrimination  of  earthquakes 
and  explosions  at  regional  distance  may  be  possible  by  using  the  coda  method  for  the 
high  frequency  case.  For  the  low  frequency  case,  we  must  consider  the  crustal  effects. 

1.  INTRODUCTION 

Numerous  studies  on  discrimination  between  underground  explosions  and 
earthquakes  at  regional  distances  during  the  last  decade  have  been  based  on  variations  of 
spectral  characteristics  of  direct  wave  phases  (Pn,  Pg,  Lg,  Rg,  etc.)  or  on  spectral  ratios 
between  different  frequencies.  Explosions  and  earthquakes  have  yielded  different 
spectral  amplitudes  for  the  same  direct  wave  phase  at  regional  distances  in  many  cases 
(e.g.,  Murphy  and  Bennett,  1982;  Pomeroy  et  al.,  1982;  Bennett  and  Murphy,  1986; 
Bennett  et  al.,  1989;  Smith,  1989)  and  can  be  used  to  discriminate  underground 
explosions  from  earthquakes.  Since  these  direct  waves  are  strongly  affected  by  a 
particular  source-receiver  path,  determining  whether  the  observed  differences  in  these 
direct  signals  are  due  to  the  source  or  propagation  path  effect  can  be  difficult,  especially 
when  the  locations  of  the  explosions  and  earthquakes  are  not  close.  Furthermore,  the  late 


P  coda  waves,  direct  S  wave  and  surface  waves  are  usually  mixed  together,  and  it  is 
difficult  to  isolate  individual  wave  phases. 

The  above  problems  can  be  solved  using  the  coda  wave  method  (Aki  and  Chouet, 
1975).  The  fundamental  characteristic  of  coda  waves  is  that  its  decay  depends  on  the 
average  properties  of  the  region  surrounding  the  source  and  station,  independent  of  any 
particular  wave  path.  This  property  allows  us  to  isolate  the  source  effect  from  path 
effects,  compare  source  differences  between  explosions  and  earthquakes,  and  develop  a 
reliable  discrimination  method. 

According  to  the  above  idea,  we  developed  a  seismic  coda  wave  discrimination 
method  and  succeeded  in  applying  it  to  discriminate  quarry  blasts  from  earthquakes  at 
local  distance  in  an  area  encompassing  the  south-central  Mojave  Desert  and  Eastern 
Transverse  Ranges  in  Southern  California.  We  found  that  the  coda  decay  rate  Qc'1  is 
significantly  higher  for  quarry  blasts  than  earthquakes  for  lower  frequencies  (f<6  Hz)  for 
lapse  time  up  to  about  30  sec.  This  result  is  attributed  to  the  greater  contribution  of 
surface  waves  to  quarry  blasts  due  to  their  shallow  source  depth.  The  difference  in  Qc_1, 
however,  disappears  for  lapse  time  greater  than  30  sec  for  the  same  frequency  range  as 
well  as  for  higher  frequencies  (f>6  Hz)  for  lapse  time  greater  than  20  sec,  suggesting  that 
the  coda  waves  at  high  frequencies,  or  at  low  frequencies  for  late  lapse  time,  are 
dominated  by  the  same  body  waves,  probably  S  waves,  for  both  quarry  blasts  and 
earthquakes.  The  coda  power  spectrum  (after  correction  of  path  attenuation)  of  the 
quarry  blasts  is  significantly  different  from  that  f  the  earthquakes.  The  curves  of  coda 
power  spectra  versus  frequency  for  quarry  blasts  decrease  more  sharply  than  earthquakes 
at  high  frequencies,  indicating  a  lack  of  energy  in  high  frequencies  for  quarry  blast 
sources  as  compared  to  earthquake  sources.  The  different  frequency  dependence  of  coda 
power  spectra  between  quarry  blasts  and  earthquakes  is  attributed  to  their  different  source 
properties  and  can  be  used  for  seismic  discrimination. 

To  extend  this  idea  to  the  seismic  discrimination  of  explosions  and  earthquakes  at 
regional  distances,  we  need  to  interpret  the  scattering  properties  of  body  waves  mixed 
with  surface  waves.  We  shall  first  develop  a  theory  for  surface  waves  scattering  and 
coupling  with  body  wave.  Using  this  theory  we  shall  investigate  the  possibility  of  the 
seismic  discrimination  of  explosions  and  earthquakes  at  regional  distance. 


2.  SEISMIC  CODA  WAVE  DISCRIMINATION  METHOD  FOR  LOCAL  EVENTS  AND  ITS 
APPLICATION  TO  QUARRY  BLAST  -EARTHQUAKE  DISCRIMINATION. 

2.1  The  Method 

Generally,  the  observed  coda  wave  decay  rate  is  very  stable,  independent  of 
source-receiver  location,  whereas  the  absolute  amplitude  of  the  coda  waves  is  dependent 
on  source  and  site  (Aki,  1969;  Aki  and  Chouet,  1975).  The  power  spectrum  of  the  coda 
wave  P(dt)  can  be  considered  as  a  product  of  three  factors: 

F(  o)l  t)=source( «))  •  site(  to)  •  path(  col  t) ,  ( 1 ) 

where  to  is  the  circular  frequency  and  t  is  the  lapse  time  measured  from  the  event  origin 
time.  Using  the  first  order  body  wave  scattering  model  of  the  coda  waves  (Aki  and 
Chouet,  1975;  Zeng  et  al.,  1990),  the  above  equation  can  be  written  explicitly  as 

lntPCcolQt^CCcoPcoQc1!,  (2) 

where  C(co)  is  determined  by  the  source  and  site  effects  at  a  given  co.  Qc* 1  is  a  parameter 
measuring  the  coda  wave  decay  rate. 

Since  the  first  order  scattering  model  neglects  the  higher  order  multiple  scattering, 
it  is  not  surprising  that  the  linear  relationship  expressed  in  Equation  (2)  holds  well  for 
only  a  finite  range  of  lapse  time.  We  observed  the  coda  decay  rates  decreasing  gradually 
with  increasing  lapse  time.  In  order  to  incorporate  the  lapse  time  dependent  coda  Qc*1, 
we  define  it  by  the  following  equation: 

~{ln[P(o)lr)r2]}  -  -<oQr\w,t),  (3) 

where  Q~\(o,t)  was  determined  by  a  moving  lapse  time  window  (centered  at  t)  in  the 
data  processing.  In  order  to  find  the  expression  for  P(colt)  by  integrating  Equation  (3),  it 
is  necessary  to  introduce  a  reference  time  tr .  As  can  be  seen  later  when  we  show  the 
observed  Qc*1  as  a  function  of  t,  we  shall  find  that  Qc*1  becomes  approximately  constant 
for  t  greater  than  a  certain  time,  which  depends  on  frequency.  It  is  then  convenient  to 
define  the  time  beyond  which  Qc_1(co,t)  becomes  constant  as  the  reference  time  tj- .  The 
integration  of  Equation  (3)  from  tr  to  t  (t^t)  gives  rise  to 


P(wl tr )t2  exp((oQ^rtr)  =  P(u)\t)t2 exp(ajQcrt)  (t*tr) 


(4) 


where  Q^,1  represents  constant  Qc'*  after  time  tr .  Equation  (4)  shows  that 

P(o/l/r  )tr2  exp(a}Q'c]tr),  t>tr,  is  a  constant  for  a  given  frequency.  The  physical  meaning 

of  the  constant  is  clear  -  it  is  the  coda  power  spectrum  corrected  for  geometrical 
spreading  and  attenuation  back  to  the  origin  time.  We  define 

Po(o,)  =  P(cohK2exp(o)Q^)  (t  *  tr)  (5) 

Obviously,  the  power  spectrum  P0((u)  we  defined  in  Equation  (5)  is  not  only 
controlled  by  the  source,  but  is  also  affected  by  the  recording  site.  If  we  compare  the 
power  spectrum  P0((o)  of  quarry  blasts  and  earthquakes  at  the  same  station,  however,  the 
difference  between  them  should  be  attributed  to  the  source  effect  only.  Integrating 
Equation  (3)  again  from  an  observation  time  t  to  the  reference  time  tr,  we  find 

ln[P(<whr)/r2]  -  ln[P((o\t)t2]  -  -cof'  QT\(o,t)dt  .  (6) 

Combining  Equations  (5)  and  (6),  we  obtain 

ln[P0(co)]  =  ln[P(«)i  t)t2]  +  ccQcrtr  -  wf'  Q t^dl  (7) 

Formula  (7)  is  obtained  under  the  assumption  that  Q ~cl{(o,f)  is  independent  of  t  for  t  a  tr . 
for  the  purpose  of  the  comparison  of  the  coda  source  power  spectra  between  explosions 
and  earthquakes,  it  can  be  generalized  to  the  case  in  which  the  coda  decay  of  explosions 
and  earthquakes  are  the  same  (instead  of  the  same  constant  Qc  after  tr). 

2.2  The  Result  of  Discrimination  Quarry  Blasts  from  Earthquakes 

A  total  of  58  earthquake  and  48  quarry  blasts  distributed  in  an  area  encompassing 
the  south  central  Mojave  Desert  and  Eastern  Transverse  Range  (Figure  1)  were  used  for 
this  study.  These  events  were  record  from  1981  to  1987  by  the  USGS-CALTECH 
Southern  California  seismograph  network.  The  quarries  are  operated  by  the  South¬ 
western  Portland  Cement  Company.  The  blasting  practice  during  this  period  was  by 
single  instantaneous  shot  in  holes  at  depths  of  about  35  feet  (J.  Stafford,  personal 
communication,  1990).  The  earthquake  focal  depths  range  from  0  to  12  km  and  their 
magnitudes  range  from  2.1  to  3.5. 


Data  from  both  quarry  blast  and  earthquake  sources  were  processed  by  an 
identical  procedure  to  study  the  possible  differences  in  coda  attenuation.  Figures  2a 
through  2d  show  the  measured  Q'1  of  earthquake  plotted  as  a  function  of  lapse  time  for 

frequencies  1.5, 3.0, 6.0,  and  12.0  Hz,  respectively,  while  Figure  2a  through  2h  show  the 
corresponding  results  for  quarry  blasts.  Each  small  open  circle  on  the  plots  represents 
one  measurement  for  a  particular  seismogram  with  a  moving  lapse  time  window  of  fixed 
width,  which  is  34  sec  for  frequency  1.5  Hz,  25  sec  for  3.0  Hz,  and  20  sec  for  6.0  and 
12.0  Hz.  The  resultant  Q~‘  values  for  all  earthquakes  were  plotted  together.  The  solid 

lines  connect  the  mean  points  (solid  circles)  calculated  by  averaging  the  individual 
measurements  in  each  8  sec  time  interval  with  4  sec  overlapping  at  the  adjacent  mean 
points.  The  standard  error  of  the  mean  is  also  shown  for  each  mean  point. 

From  the  comparison  of  Figures  2a  through  2h,  we  find  a  significant  difference  in 
Q”1  between  quarry  blasts  and  earthquakes  for  the  lower  frequencies  of  1.5  and  3  Hz  for 
lapse  times  less  than  30  sec.  this  difference  in  Q’1 ,  however,  disappears  for  lapse  time 
greater  than  30  sec  in  the  same  frequency  range  as  well  as  for  higher  frequencies  6  and  12 
Hz  for  lapse  time  greater  than  20  sec. 

We  interpret  the  above  observation  as  follows:  since  quarry  blasts  are  shot  at 
extremely  shallow  depth  (35  feet),  they  generate  more  surface  waves  than  earthquakes. 
These  surface  waves  attenuate  quickly  because  their  propagation  path  is  confined  to 
shallow  depths  where  seismic  attenuation  is  strong.  Body  waves  thus  become  dominant 
over  surface  waves  for  t>20  sec  at  6  and  12  Hz  and  for  t>30  sec  at  1.5  and  3.0  Hz.  Once 
body  waves  dominate  over  surface  waves,  the  coda  decay  rate  becomes  the  same  for 
quarry  blasts  and  earthquakes. 

The  above  observed  coda  decay  rate  between  quarry  blasts  and  earthquakes  leads 
us  to  our  basic  idea  for  discrimination,  namely,  the  comparison  of  their  power  spectra 
P0(cd)  inferred  at  a  time  later  than  tr  at  a  single  station.  These  power  spectra  P0(co)  will 
not  be  influenced  by  surface  waves  and  will  not  depend  on  media  properties,  since  body 
waves  dominate  the  late  coda  waves  of  both  quarry  blasts  and  earthquakes,  with  the  same 
scattering  and  attenuation  properties.  We  have  chosen  stations  SIL  and  SDW  for  our 
power  spectra  comparison,  because  these  two  stations  were  operated  satisfactorily 
through  the  time  period  from  which  our  data  were  collected  and  were  located  close  to 
both  the  quarry  blasts  and  earthquakes  used  in  this  study  (see  Figure  1).  Figure  3  shows 


the  coda  power  spectrum  P0(to)  vs.  frequency  for  station  SIL  for  earthquakes  and  quarry 
blasts,  respectively.  Each  curve  in  these  plots  represents  one  event.  We  only  plot  those 
curves  for  which  P0(o>)  is  available  for  at  least  three  frequency  bands.  All  curves  were 

normalized  to  the  same  value  at  frequency  1.5  Hz  for  convenience  of  comparison.  In 
order  to  test  the  assumption  of  constant  Q"1  for  t  s  tr ,  we  also  used  a  different  reference 

time  tr  with  tr=20, 32, 40, 52  sec  corresponding  to  Figures  3a,  3b,  3c,  and  3d  for 
earthquakes  and  Figures  3e,  3f,  3g,  and  3h  for  quarry  blasts.  We  can  see  that  P0(o))  in 
Figure  3b  is  very  similar  to  those  in  figures  3c  and  3d.  The  same  feature  is  also  observed 
in  figures  3f,  3g,  and  3h  for  quarry  blasts.  The  shape  of  the  power  spectrum  P0(cd)  is  not 
affected  very  much  by  the  choice  of  tr  as  long  as  tr>30  sec,  which  confirms  the  validity  of 
our  assumption  that  Q'1  is  constant  for  t  s  tr .  By  comparing  earthquake  power  spectrum 
P0(to)  with  that  of  quarry  blasts  at  the  same  reference  time  tr  (such  as  Figures  3a  and  3e, 
3b  and  3f,  etc.),  we  see  that  there  is  a  significant  difference  between  these  two  groups. 

The  power  spectra  P0(o>)  of  quarry  blasts  are  decreasing  more  quickly  than  earthquakes 
with  increasing  frequency,  indicating  that  the  quarry  blast  source  generates  less  high 
frequency  energy  than  earthquakes.  The  same  comparison  was  also  made  between 
earthquakes  and  quarry  blasts  at  station  SDW.  The  results  are  consistent  with  those  from 
station  SIL.  This  difference  in  source  spectra  frequency  dependence,  therefore,  can  be 
used  for  seismic  discrimination  between  quarry  blasts  and  earthquakes. 

To  check  the  depth  and  magnitude  effects  on  the  power  spectrum  of  earthquakes 
and  quarry  blasts,  we  also  chose  one  of  the  shallower  earthquakes  (focal  depth  of  0.6  km), 
one  of  the  deeper  earthquakes  (focal  depth  of  1 1.9  km),  and  one  of  the  quarry  blasts  for 
comparison.  The  magnitude  of  these  three  events  are  nearly  the  same  (2.4  ~  2.5).  The 
seismograms  of  the  quarry  blast,  the  shallow  earthquake  and  the  deeper  one  for  station 
SIL  are  shown  in  Figures  4a,  4b  and  4c.  The  top  traces  in  these  Figures  show  the 
unfiltered  seismograms  while  the  other  traces  show  the  results  of  band  pass  filtering  of 
the  seismograms  with  center  frequencies  of  1.5, 3,  6,  and  12  Hz.  The  maximum 
amplitude  of  these  three  events  are  almost  the  same  from  a  comparison  of  their  unfiltered 
seismograms.  However,  for  the  band  pass  seismograms  at  12  Hz,  the  quarry  blast 
amplitude  is  significantly  smaller  when  compared  to  the  two  earthquakes. 

The  time  dependent  spectral  characteristics  of  these  three  events  are  plotted  in 
Figures  5a,  5b  and  5c.  The  contour  lines  in  these  figures  show  the  iso-amplitude 
spectrum  line  plotted  in  the  frequency-time  diagram.  The  frequency  range  is  1  to  25  Hz, 


and  the  lapse  time  range  is  from  0  to  50  sec.  The  frequency  -time  dependence  of  the 
amplitude  spectra  contour  is  distinctly  different  between  the  earthquakes  and  the  quarry 
blast.  The  spectra  for  the  quarry  blast  shows  very  rapid  decay  with  increasing  frequency 
while  those  for  earthquakes  cover  much  broader  spectral  ranges.  Figure  6  shows  the 
power  spectrum  P0(oo)  vs.  frequency  for  the  deep  earthquake  (line  1),  the  shallow  one 
(line  2)  and  the  quarry  blast  (line  3).  The  depth  difference  between  the  quarry  blast  and 
the  shallow  earthquake  is  about  a  half  km  while  the  difference  between  the  shallow 
earthquake  and  the  deep  one  is  about  1 1  km.  Thus,  we  can  conclude  that  the  depth  effect 
is  not  the  main  reason  for  causing  the  source  spectrum  shape  difference  between  quarry 
blasts  and  earthquakes  in  our  results. 

Therefore,  we  can  conclude  that  the  observed  source  spectrum  difference  between 
earthquakes  and  quarry  blasts  are  mainly  due  to  different  source  processed,  such  as  the 
difference  in  source  mechanism  and  the  interaction  between  source  and  near-source 
media.  This  is  physically  reasonable  because  the  quarry  blasts  are  shot  at  shallow  depths 
and  the  near-source  media  are  less  consolidated  than  the  media  generating  earthquakes, 
which  may  result  in  a  lack  of  high  frequency  energy  generated  compared  to  the 
earthquakes. 

3.  Seismic  Scattering  Theory  for  surface  waves  and  Body  waves 

AND  ITS  IMPLICATION  OF  DISCRIMINATION  TO  EXPLOSION  -E  ARTHQUAKE  AT 

REGIONAL  DISTANCE 

In  the  preceding  section,  we  have  shown  that  at  local  distances,  if  we  can  properly 
remove  the  path  effects  from  the  observed  power  specter,  we  can  discriminate  the  quarry 
blast  and  the  earthquake  from  the  shape  of  the  source  power  spectra,  for  local  events, 
coda  Qc1  method  provides  an  effective  mean  to  isolate  the  source  effect  from  path 
effects  as  shown  above.  To  extend  this  discrimination  method  to  the  discrimination  of 
explosions  and  earthquakes  at  regional  distances,  we  need  a  proper  means  to  remove  the 
path  effects  at  regional  distances.  Since  shallow  events  produce  significant  surface 
waves,  a  complete  scattering  theory  must  include  the  scattering  processes  of  surface 
waves,  and  the  coupling  process  between  body  waves  and  surface  waves.  In  this  section, 
we  shall  develop  an  energy  transfer  theory  for  a  mixture  of  surface  wave  and  body  wave 
to  study  their  scattering  and  coupling  processes.  Using  this  theory  we  shall,  first,  explain 
the  physical  basis  for  successful  discrimination  of  quarry  blasts  and  earthquakes  at  local 
distances  described  in  the  preceding  chapter.  Then,  by  fitting  our  synthetic  coda 


attenuation  curves  with  the  observed  data,  we  will  estimate  the  scattering  and  coupling 
coefficients.  Finally,  with  the  derived  scattering  and  coupling  coefficients,  we  shall 
simulate  the  surface  wave  and  body  waves  scattering  processes  at  regional  distances  and 
investigate  the  feasibility  of  discrimination  of  explosions  and  earthquakes  at  regional 
distance.by  using  the  coda  Qc_1  method. 

3.1  Energy  transfer  theory  for  Rayleigh  wave 

In  this  section,  we  shall  consider  a  simple  case  in  which  the  background  medium 
is  a  homogeneous  half-space.  The  only  surface  wave  in  this  case  is  the  Rayleigh  wave. 
For  simplicity,  we  neglect  the  conversions  between  the  body  waves  and  the  Rayleigh 
wave.  Then  the  energy  distribution  of  Rayleigh  wave  in  an  absorption  and  scattering 
medium  at  an  arbitrary  point  x  for  frequency  o>  can  be  written  as, 

E^x.co)  =  G(x,xs)exp[-qV-rsl]e§(xs,to)  +  J  cKXj)G(x,x  Jexp[-qV-rjllER(Xj«>), 

j 

(8) 

Where,  r=  x-ez(x*ez)  is  a  horizontal  vector,  and  qR  is  the  attenuation  coefficient  that 
includes  the  intrinsic  (qf )  and  scattering  attenuation  (qf ),  i.e.,  qR  =  qf  +  qf.  In  eq.(8), 
E^x.co)  is  the  seismic  energy  per  unit  volume  carried  by  Rayleigh  wave  at  point  x  and 
frequency  m.  The  first  term  in  the  right-hand-side  of  (8)  represents  the  direct  Rayleigh 
wave  energy  radiated  from  source,  and  the  second  term  in  the  right-hand-side  of  (8) 
describes  the  total  scattering  energy  from  all  scatters  (xj).  Where  a  is  the  scattering 
cross-section,  G(x,Xj)  and  expf-q^r-rjl]  describes  the  geometrical  spreading  and 
attenuation,  respectively.  It  should  be  noted  that  the  attenuation  is  caused  only  by  the 
horizontal  propagation,  since  the  surface  wave  only  propagates  horizontally.  For  weak 
scattering  case,  the  scattering  cross-section  a  can  be  analytically  derived  by  using  Bom 
approximation  (Aki  and  Richards,  1980).  Here,  for  simplicity,  we  assume  that  a  is  a  con¬ 
stant.  The  geometrical  spreading  function  G(x,x'),  however,  can  be  determined  by  the 
consideration  of  energy  conservation. 

Geometrical  Spreading  Function  of  Rayleigh  Wave: 

The  direct  Rayleigh  wave  (Rayleigh  wave  in  the  background  medium)  has  the 
following  form  (Aki  and  Richards,  1980), 


(9) 


UR(x,<o)  -  Pr(X/)[Pr(Xs),F(Xs)]  exp{i [kRlr-rs I+ji/4] }, 

V  (jt/2)kRlr-rsl 


where. 


PR(x)  =  rx(z,cD)er  +  ir2(z,(o)ez , 


with  er  =  (r-rs)/lr-rsl.  The  corresponding  energy  density  is 

E*(x,o>)  =  WlU^x.o))!2  =  P^2.  -R.(^li|PR(xs)»F(xs)l2  .  (10) 

2  JtkR  |r_rs| 

For  steady  state,  the  energy  rate  should  be  a  constant  for  any  closed-surface  that  encloses 
the  source,  namely. 


|  (ES(x, 


a))vR}*? dS  =  constant,  for  xs  inside  S. 


Where,  s  is  the  normal  vector  of  the  surface  element  dS,  and  vR  is  the  velocity  vector  of 
the  Rayleigh  wave.  Since  the  surface  wave  propagates  along  the  horizontal  direction,  we 
consider  a  cylindrical  surface  that  leads  to  a  simple  surface  integration. 


{E£(X,0))Vr} 


Jp+oo  [2 

dz  I 
0  Jo 


Rd0  ES(x,cd)vr 


Where,  R  =  lr-rsl.  Substituting  eq.(10)  in  (11),  we  find 


Eq(xs)  =i-j>  {E£(x,cd)vr}-s  dS  =  IPR(z)l2dz  IPR(xs)»F(xs)l2 

2  f” 

=  IF(xs)l2{lri(Zs)nrl2  +  lr2(Zs)nzl2}  {lri(z)l2  +  lr2(z)l2}  dz 

kR  Jo 


Using  the  above  result,  eq.  (10)  can  be  rewritten  as 


Eo(x,cd) 


gR(z,co) 

2jilr-rsl 


«&(*»)> 


(13) 


where. 


gR(Z,«>) 


IPr(z)I2 


IPR(z')l2dz' 


[ri(z,to)2  +  r2(z,o>)2] 
[ri(z',o>)2  +  r2(z',a>)2]dz' 


(13a) 


Thus,  we  obtain  the  geometrical  spreading  function  of  Rayleigh  wave  energy  propagation 
as 

G(x,x’)  =  i^^.  (14) 

2jtlr-r'l 


From  equation(6a),  we  can  verify  the  following  identity. 


gR(z,(o)dz  =  1 . 


(14a) 


Integral  equation  of  scattering  Rayleigh  wave  in  random  medium: 

We  are  now  ready  to  set  up  the  basic  integral  equation  for  Rayleigh  wave 
scattering  in  random  medium.  Substituting  eq.(14)  in  (8),  and  assuming  that  the  random 
scatters  are  uniformly  distributed  and  can  be  described  by  a  continuous  distribution 
tio(density  of  scatters),  we  obtain  the  following  integral  equation, 

ER(x,ca)  =  gR(Z’t°)e"TlR|r:r-eg(x,(o)  +  j  r]R  M^e-^r-r'iER(x,a))dV(x') ,  (15) 

2jtlr-rsl  I  2jilr-r'l 

where,  =  iioO.  Eq.(15)  is  similar  to  the  energy  transfer  equation  for  body  wave  (Wu, 
1985;  Zeng  et  al.,  1991),  but  with  different  geometrical  spreading  functions  and  the 
attenuation  factors  that  are  caused  by  the  characterizations  of  surface  wave  propagation. 
The  dependence  of  lr-r'l  indicates  the  horizontal  propagation  property  of  the  surface 


wave,  whereas,  the  gR(z,co)  explores  the  depth  distribution  of  surface  wave  energy. 
Integral  equation  (15)  can  be  solved  by  using  spatial  domain  Fourier  transform  method. 
First,  we  note  that  equation  (15)  can  be  rewritten  as 


ER(x,o>)  =  gR(z,co)E2(r,o)). 


(16) 


Inserting  equation  (16)  into  (15)  and  using  the  identity  (14a),  we  obtain 

e-r]Rlr-rJ 


Ef(r,o))  = 


2jtlr-r, 


■£o(Xs,(0)  + 


I  T,Re-^r-f'1 
Js  s  2jtlr-rl 


E?(r,a>)  dZ(r). 


(17) 


Here,  2(r)  =  {(x,y)l-o°<x<+°°,  -x<y<+oo}.  The  corresponding  Fourier  transform  is 

B(k,a>)  =  G2(k)^(xs,co)  +  T£05(k)g(k,a>).  (18) 

Where,  we  assumed  rs  =  0,  and  G2(k,<o)  is  given  by 


+o°  ^ 

G^(k,o))=  I  (  dxdy  £^-exp(ik*r)  =  I  J0(kr)e-TiRrdr  =  — — L___ 

I)  2ltr  Jo  Vk2  +  (rf)2 


(18a) 


Substituting  this  result  in  eq.(18),  then  taking  the  inverse  Fourier  transform  over  k,  we 
finally  obtain 

ER(x,oa)  =  gR(z,o))PR(r,rf  ,qf)  j£(xs,o )) ,  (19) 

where  PR(r,qR,rif)  is  the  inverse  Fourier  transform  of  G^(k,o>), 


PrO\t|r,tir) 


Jo(kr) 


kdk 


Vk2  +  (qR)2 


(19a) 


Solution  (19)  has  a  clear  physical  meaning.  The  term  eR(xs,(o),  as  defined  earlier, 
contains  the  seismic  source  information.  Function  gR(z,co)  represents  the  depth 
dependence  of  Rayleigh  wave  energy,  whereas  the  function  PR(r,qR,qR)  describes  the 
propagation  and  attenuation  processes  of  the  Rayleigh  wave  energy.  PR(r,T]R,T|f)  can  be 
obtained  by  evaluating  integral  (19a).  For  a  pure  absorption  medium  (r^O),  we  find 


PR(r,qR,0)  =  e-’A/Zjir .  This  is  consistent  with  our  direct  Rayleigh  wave  (solution  in 
background  medium).  For  general  absorption  and  scattering  media,  we  can  numerically 
evaluate  the  propagation  and  attenuation  function  PR(r,qR,qR). 

It  is  noted  that  equation  (15)  describes  a  stationary  energy  transfer  process,  i.e., 
Er(x,co)  represents  the  amplitude  of  the  spectrum  of  scattered  Rayleigh  wave.  To  obtain 
the  time-history  of  the  scattered  Rayleigh  energy  for  a  given  frequency  go  and  the 
observational  point  x,  we  introduce  a  time-delay  phase  factor  of  exp[-/Qlr  -  r'  I  lcR ] ,  then 
take  inverse  Fourier  transform  over  frequency  Q.  During  the  inverse  Fourier 
transformation,  the  frequency  go  is  kept  as  a  constant  which  is  defined  as  center 
frequency.  Thus  we  can  obtain  the  time-history  of  scattered  energy  for  a  given 
frequency  go.  Figure  7  shows  the  energy  distributions  for  the  case  of  qf  =0.01  km-1, 
rjf  =0.02  km1  and  co=l  H£,.  As  we  expected  the  scattered  energy  arrived  at  time  of 
t=lr-r'l  lcR. 

3.2  Energy  transfer  theory  for  S  and  Rayleigh  waves’  scattering  processes 

In  the  preceding  section,  we  have  considered  only  the  scattering  processes  of 
surface  waves,  and  neglected  the  body  waves’  conversions.  In  this  section,  we  shall 
consider  the  coupling  effect  of  body  waves  with  surface  waves.  For  simplicity,  we 
consider  only  the  coupling  process  between  the  S  wave  and  Rayleigh  wave.  The  energy 
transfer  equation  described  such  an  coupled  scattering  process  can  be  written  as  follows, 

ER (x, co)  =  Gr( x,xs)  ■  e~ 'eR (xs,w) 

+fvGR(x,x')-e-v**-''{r1™ER(x\(o)  +  ri°REs(x',(0)}dV(ji!),  (20a) 

Es(x,a>)  =  Gs(x,xs)  •e-”S|r-r'i£os(xJ,co) 

+JvGs(x,x')  ER(x',w)  +  tj*sEs(x!  ,eo)}rfV(x')  ,  (20b) 

where, 

Es(x,oo):  seismic  energy  carried  by  S  wave  per  unit  volume  at  x ; 

E^x,  w):  seismic  energy  carried  by  Rayleigh  wave  per  unit  volume  at  x; 

r]R  :  total  attenuation  coefficient  for  Rayleigh  wave,  qR  =  r]R  +  +qfR; 

r|!p:  Rayleigh  to  Rayleigh  waves  scattering  coefficient; 

qf R :  S  to  Rayleigh  waves  scattering  coefficient; 
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tjs  :  total  attenuation  coefficient  for  S  wave,  r\s  =  rjp  +  +rifs; 

r|f  :  absorption  coefficient  for  S  wave; 

riP :  Rayleigh  to  S  waves  scattering  coefficient; 

r|fs  :  S  to  S  waves  scattering  coefficient; 

8o(xs,co):  rate  of  S  wave  energy  rate  radiated  from  source  divided  by  the  velocity 
of  S  wave; 

£*(xs,(o ) :  rate  of  Rayleigh  wave  energy  radiated  from  source  divided  by  velocity 
of  Rayleigh  wave; 


Gs(x,x'):  the  geometrical  spreading  function  for  S  wave  in  a  half-space  medium, 

and  it  can  be  approximated  by 

Gs(x,x')  =  -  l  -  -■■+ - 1 — -  , 

4rclx-x  I  43tlx-(x')  I 


where,  (x’)*  is  the  image  point  of  x’  with  respect  to  the  free  surface  z=0,  i.e.,  (x’)*=(x\ 

y’.-z’)- 

Equations  (20a)  and  (20b)  can  be  further  simplified.  Making  a  two-dimensional 
Fourier  transform  over  the  horizontal  variable  (x,  y),  equation  (20a,  b)  can  be  reduced  to, 

E?(K,z,o))  =  GR(K,z;zs;o))£o(Xg,ft))  (21a) 


Jm+  co 

GR(K,z;z';co){riP^(K,z»  +  r,f  ES(K,z',co)}dz', 

0 

Ef(K,z,co)  =  Gs(K,z;zs;(o)eS(xs,(o) 

JA+CO 

Gs(K,z;z';co){tisSEr(K,z;z,;(o)  +  Tif^(K,z',(o)}dz' , 

0 

where,  K  is  the  horizontal  wave-number  and 
GS(K, 


(21b) 


:,z;z';o))=  I 

Jo 


and 


G(r,z;z';to)exp(-T]sVr2+(z-z')2  )J0(Kr)rdr, 


Solving  the  coupled  equations  (21a)  and  (21b),  we  can  obtain  the  solution  of  energy 
transfer  processes  of  S  and  Rayleigh  waves.  The  time-history  of  total  scattering  energy 
can  be  determined  by  introducing  proper  time-delay  phase  factors  to  these  coupled 
equations.  Figure.  (8)  show  the  time  history  of  total  scattering  energy  of  S  and  Rayleigh 
waves,  for  the  case  of  T|f=0.01(km)_1,  qP=0.01(km)‘1,  qp=0.02(km)"1,  rjf^O.OSikm)1, 
r^RR=0.02(km)-1  and  O.Ol(km)-1  for  various  source  depths.  Our  results  indicate  that 
the  shallower  the  source  depth  locates  or  the  lower  the  frequency  is,  the  contribution  of 
scattering  Rayleigh  waves  is  larger.  This  is  consistent  with  the  observed  results  of  Su  and 
Aki  (1991)  as  shown  in  Figure  2,  where  earthquakes  have  an  averaged  focal  depth  of  8 
km,  whereas  the  averaged  focal  depth  of  quarry  blasts  is  about  10m. 

To  quantitatively  interpret  the  observation  shown  in  Figure  2,  and  determine  the 
coupling  coefficients,  we  shall  compare  our  synthetic  attenuation  curves  with  the 
observed  ones.  The  synthetic  attenuation  curves  are  calculated  using  the  same  method  as 
described  in  section  2. 1  from  our  theoretical  scattered  S  wave  energy  envelope.  We  here 
consider  only  the  S  wave  coda  attenuation,  Qc'* ,  that  mainly  includes  the  S  to  S  waves 

scattering  as  well  as  Rayleigh  to  S  waves  scattering.  For  simplicity,  we  neglect  the  S  to 
Rayleigh  conversions,  thus  the  only  coupling  coefficient  needs  to  be  determined  is  rj™. 

Physically,  the  discrimination  between  the  earthquakes  and  quarry  blasts  can  be  attributed 
to  the  contributions  of  Rayleigh  wave,  namely,  the  measure  of  coda  attenuation  of  quarry 
blasts  contains  significant  Rayleigh  wave’s  conversions  because  of  the  shallow  focal 
depth,  while  the  coda  attenuation  of  the  earthquakes  contains  little  Rayleigh  waves 
conversions  because  of  the  deeper  focal  depth.  Based  on  this  point  we  shall,  first, 
determine  the  S  wave’s  intrinsic  and  scattering  coefficients,  i.e.,  rjf  and  rj*s,  by  matching 

our  synthetic  S  wave’s  attenuation  with  the  observed  earthquake’s  coda  attenuation,  Qc_1. 
At  this  stage  we  let  %s=0.  Having  the  coefficients  of  rjf  and  rjsss,  we  then  determine  the 
coefficients  of  r/f5,  rj™  and  rjf  by  fitting  our  synthetic  attenuation  curves  with  the 
observed  quarry  blasts  coda  attenuation,  Qc'f 

After  more  than  20  trials,  we  found  the  attenuation  coefficients  for  the  best 
matching  with  the  observed  data  as  given  in  Table  1.  And  the  corresponding  coda  atten¬ 
uation  curves  are  shown  in  Figure  9.  Our  results  indicate  that  the  coupled  scattering 


coefficients,  rjf,  which  describes  the  scattering  conversion  from  Rayleigh  wave  to  S 

wave,  are  strongly  depends  on  frequency  in  the  following  way:  for  frequency  lower  than 
3  Hz,  rj™ =0.05  km-1;  while  for  frequency  greater  than  6  Hz  r/f  <0.0005  km'1.  Therefore, 

there  may  exist  a  critical  frequency  between  3  Hz  and  6  Hz  beyond  which  the  Rayleigh 
wave  conversions  are  negligible,  while  below  which  the  Rayleigh  wave  conversions  to  S 
wave  are  significant. 

3.3  Discussion  on  the  Discrimination  of  Earthquakes  and  Explosions  at  Regional 

Distance 

Having  the  theoretical  studies  on  the  discrimination  of  earthquakes  and  quarry 
blasts  at  local  distance,  we  can  discuss  the  possibility  of  the  discrimination  of  earthquakes 
and  explosions  at  regional  distance  using  the  coda  method.  In  the  following  discussion 
we  assume  that  the  intrinsic  and  scattering  attenuation  coefficients  of  S  and  Rayleigh 
waves,  and  the  coupled  scattering  coefficients  obtained  in  the  above  study  are  still  valid 
for  the  events  at  regional  distance.  Using  those  derived  coefficients  we  calculated  the 
synthetic  coda  attenuation  curves  for  earthquakes  and  explosions  at  regional  distances  as 
shown  in  Figure  10.  Here,  we  assume  that  the  hypocenter  is  600  km,  and  the  focal  depth 
is  10  km  for  earthquakes  and  0.01  km  for  explosions.  As  shown  in  Figure  10,  our 
theoretical  simulations  show  that  the  Rayleigh  wave  conversion  to  scattered  S  wave  at 
regional  distance  is  much  more  significant  than  that  at  local  distance,  because  the 
Rayleigh  wave  has  a  weaker  geometrical  spreading  factor  that  S  wave.  For  lower 
frequency  (e.g.,  f=3  Hz  or  1.5Hz),  the  synthetic  result  predicts  a  strong  increase  in  coda 
attenuation  that  starts  just  after  the  Rayleigh  wave's  arrival  and  last  about  50  seconds. 

For  higher  frequency  (e.g.,  f=6  Hz  or  12  Hz),  a  similar  increase  is  also  predicted  by  our 
synthetic  result,  although  it  is  not  as  strong  as  the  one  for  lower  frequencies.  It  should  be 
pointed  out  that  this  increase  in  the  coda  attenuation  estimate  for  lower  frequency  is  an 
upper-bound,  the  actual  predicted  increase  in  coda  attenuation  should  be  below  this 
amount;  because  the  value  of  coupled  scattering  coefficient  rj™  is  given  in  the  range  of  0 

to  0.0005  km'1. 

In  order  to  compare  the  above  result  with  the  observation,  we  have  analyzed  the 
regional  events  recorded  by  the  Chinese  Digital  Seismic  Network  (CDSN).  Figure  1 1 
shows  the  seismograms  of  the  earthquakes  and  explosions  used  in  this  analysis.  These 
seismograms  are  recorded  by  the  CDSN  at  AKK  station.  The  results  of  our  analysis,  i.e., 
the  coda  Qc  1  of  S-wave  for  both  the  explosions  and  earthquakes  are  plotted  in  Figure  12, 


which  show  a  similar  result  as  predicted  only  for  the  high  frequency  case.  The  reason 
why  our  theoretical  predictions  fail  for  the  low  frequency  case  is  due  to  the  absence  of  the 
crustal  layer  in  our  model.  When  frequency  becomes  low,  the  corresponding  wavelength 
becomes  comparable  to  the  thickness  of  the  crust,  therefore,  the  crustal  effects  on  the 
seismic  scattering  cannot  be  neglected.  In  the  future  study,  we  shall  consider  the  crustal 
effects  on  the  seismic  waves  scattering  processes,  in  order  to  obtain  a  more  complete 
scattering  theory  at  regional  distances. 

4.  Discussions  and  Conclusions 

We  have  successfully  discriminated  the  quarry  blasts  and  earthquakes  at  local 
distance  from  the  observed  source  power  spectra  by  using  the  coda  Q  method.  To  extend 
this  discrimination  method  to  the  events  at  regional  distances,  we  have  developed  the 
energy  transfer  theory  for  a  mixture  of  Rayleigh  and  S  wave’s  propagation  in  a  random 
scattering  and  attenuation  medium.  Applying  this  theory  to  the  quarry  blast-earthquake 
problem  at  local  distance,  we  can  explain  the  observations  satisfactorily.  Our  study 
indicates  that  the  Rayleigh  waves'  scattering  become  important  for  the  events  with 
shallower  focal  depth  such  as  the  quarry  blasts  and  nuclear  explosions.  By  fitting  our 
theoretical  attenuation  curves  to  the  observed  ones  at  local  distances,  we  have  obtained 
the  S  and  Rayleigh  waves'  intrinsic  and  scattering  coefficients,  and  the  coupled  scattering 
coefficient  of  Rayleigh  with  S  waves.  Assuming  that  these  coefficients  are  valid  for  the 
events  at  regional  distances,  we  have  compared  the  predicted  and  observed  coda  envelope 
at  regional  distances  for  the  purpose  of  discrimination  between  the  earthquakes  and 
explosions  at  regional  distance.  Our  prediction  shows  that  the  Rayleigh  wave's 
contributions  become  more  significant  for  shallower  events  (explosions),  while  for  deeper 
events  (earthquakes)  the  Rayleigh  wave's  contributions  are  still  negligible.  Comparison 
between  predicted  and  the  data  of  CDSN  suggests  that  the  discrimination  of  earthquakes 
and  explosions  at  regional  distance  may  be  possible  by  using  the  coda  method  for  the 
high  frequency  case.  For  the  low  frequency  case,  we  must  consider  the  crustal  effects.  In 
other  words,  a  scattering  theory  for  Lg  waves  is  needed.  We  plan  to  do  this  work  in  the 
future. 
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Figured  Coda  Qc'1  vs.  lapse  time  obtained  using  the  events  shown  in  Figure  1  with  (a)  fl>)  (c)  and  fd* 
for  earthquakes  and  (e),  (f),  (g),  and  (h)  for  quarry  blasts.  Each  open  circle™  the  plot  remeants'one  ^ 
measurement  for  a  pwticular  seismogram  on  a  time  window  of  34  sec  for  frequency  1  5lKid  25  sec  for 

aJSriUff  2L  hSry6  3fd  12  H2*  ThesoiidUne  connects  the  mean  points  (solidckcles^culated  by 
averaging  the  individual  measurements  in  each  8  sec  time  interval  with  4  sec  overlapping  at  the  adjacent 
mean  points.  The  standard  error  of  the  mean  is  also  shown  for  each  mean  pnin* 
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Figure  3:  Power  spectrum  P0(<o)  vs.  frequency  obtained  at  station  SIL  with  (a),  (b),  (c),  and  (d)  correspond 
to  reference  times  of 20, 32, 40, 52  sec  for  earthquakes  and  (e),  (0,  (g)  and  (h)  for  quarry  blasts.  Each 
curve  on  these  plots  represents  one  event.  We  only  plotted  those  curves  for  which  P0((o)  is  a  vail  Ate  far  at 
least  three  frequency  lands.  All  curves  were  normalized  to  the  same  values  at  frequency  1.5  Hz  for 
convenience  of  comparison. 
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Figure  4:  MttNfams  of  (a):  a  quarry  blast  (2/13/85,  N34°  37.80’,  W1 17°  720’,  depth  0  km,  magnitude 
2.5,  epicentral  distance  41  km),  (b):  a  shallow  earthquake  (5/6/86,  N34°  20.31’,  W1 17°  5.28’,  depth  0.6 
km,  magnitude  2.4,  epicentral  distance  24  km),  (c):  a  deep  earthquake  (5/31/86,  N34°  6.43’,  W1 16°  36.58’, 
depth  1 1.9  km,  magnitude  2.4,  epicentral  distance  33  km),  recorded  at  station  SIL.  The  top  unfiltered 
seismogram  is  followed  by  bandpass-filtered  seismograms  with  center  frequencies  of  1.5, 3.6  and  12  Hz, 
respectively. 


Figure  Sc:  Contours  of  amplitude  spectra  plotted  in  a  frequency  vs.  time  diagram  for  the  deep  p^rthqimirr 
shown  in  Figure  4c. 
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Figure  6:  Power  spectrum  Po(co)  vs.  frequency  for  the  deep  earthquake  (line  1),  the  shallow  one  (line  2), 
and  the  (parry  blast  (line  3),  obtained  at  stations  SIL  at  q=52sec.  The  three  curves  were  normalized  to  the 
same  value  a  frequency  1 JS  Hz  for  convenient  of  comparison. 
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Figure  8:  Focal  depth  and  frequency  co  dependence  of  total  scattered  energy  of  Rayleigh 
and  S  wave*,  where  77*  =0.03  (km)-l,  77™  =0.02  (km)-l,  =0.01(km)*l,  T]f=0.0l(kmy\ 

f  =0.02(h«i)r1,  rj^sO.OlOan)*1  and  r=30km.  For  the  cases  of  zs  =0.0  lkm  (solid  line), 
=lkm  (dotted  line)  and  zs=10km  (dashed  line). 
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Figure  11:  Seismograms  of  the  earthquakes  and  explosions  recorded  by  CDSN  at  AKK 
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Kgure  12  shows  the  apparent  O'1  after  S-arrival,  for  earthquakes  and  explosions  at 
regional  distances,  for  different  frequencies.  ana  explosions  at 


